home *** CD-ROM | disk | FTP | other *** search
/ Complete Linux / Complete Linux.iso / docs / apps / circuits / spice2g6.z / spice2g6 / spice / Fortran / moseq3.f < prev    next >
Encoding:
Text File  |  1989-02-03  |  8.6 KB  |  334 lines
  1.       subroutine moseq3(vds,vbs,vgs,gm,gds,gmbs,
  2.      1   qg,qc,qb,cggb,cgdb,cgsb,cbgb,cbdb,cbsb)
  3.       implicit double precision (a-h,o-z)
  4. c
  5. c     this routine evaluates the drain current, its derivatives and
  6. c     the charges associated with the gate, channel and bulk
  7. c     for mosfets based on semi-empirical equations
  8. c
  9. c spice version 2g.6  sccsid=mosarg 3/15/83
  10.       common /mosarg/ vto,beta,gamma,phi,phib,cox,xnsub,xnfs,xd,xj,xld,
  11.      1   xlamda,uo,uexp,vbp,utra,vmax,xneff,xl,xw,vbi,von,vdsat,qspof,
  12.      2   beta0,beta1,cdrain,xqco,xqc,fnarrw,fshort,lev
  13. c spice version 2g.6  sccsid=status 3/15/83
  14.       common /status/ omega,time,delta,delold(7),ag(7),vt,xni,egfet,
  15.      1   xmu,sfactr,mode,modedc,icalc,initf,method,iord,maxord,noncon,
  16.      2   iterno,itemno,nosolv,modac,ipiv,ivmflg,ipostp,iscrch,iofile
  17. c spice version 2g.6  sccsid=knstnt 3/15/83
  18.       common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
  19.      1   gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox,
  20.      2   pivtol,pivrel
  21. c
  22.       equivalence (xlamda,alpha),(vbp,theta),(uexp,eta),(utra,xkappa)
  23.       data coeff0/0.0631353d0/,coeff1/0.8013292d0/,coeff2/-0.01110777d0/
  24. c
  25. c     icharg=1 causes charges to be computed
  26. c     icharg=0 bypasses the computation of charges
  27. c
  28. c     icharg=1
  29. c     if (mode.ne.1) go to 10
  30. c     icharg=0
  31. c     if (modedc.eq.2.and.nosolv.ne.0) icharg=1
  32. c     if (initf.eq.4) icharg=1
  33. c
  34. c     reference cdrain equations to source and
  35. c     charge equations to bulk
  36. c
  37. 10    continue
  38.       icharg=0
  39.       vgb=vgs-vbs
  40.       vfb=vbi-phi
  41.       vdsat=0.0d0
  42.       qg=0.0d0
  43.       qb=0.0d0
  44.       qc=0.0d0
  45.       cgdb=0.0d0
  46.       cbdb=0.0d0
  47.       onxl=1.0d0/xl
  48.       eta=eta/(xl*xl*xl)
  49. c
  50. c.....square root term
  51. c
  52.       if ( vbs.gt.0.0d0 ) go to 120
  53.       phibs=phi-vbs
  54.       sqphbs=dsqrt(phibs)
  55.       dsqdvb=-0.5d0/sqphbs
  56.       go to 200
  57. 120   continue
  58.       sqphis=dsqrt(phi)
  59.       sqphs3=phi*sqphis
  60.       sqphbs=sqphis/(1.0d0+vbs/(phi+phi))
  61.       phibs=sqphbs*sqphbs
  62.       dsqdvb=-phibs/(sqphs3+sqphs3)
  63. c
  64. c.....short channel effect factor
  65. c
  66. 200   continue
  67.       if ( (xj.eq.0.0d0).or.(xd.eq.0.0d0) ) go to 210
  68.       wps=xd*sqphbs
  69.       onxj=1.0d0/xj
  70.       xjonxl=xj*onxl
  71.       djonxj=xld*onxj
  72.       wponxj=wps*onxj
  73.       wconxj=coeff0+coeff1*wponxj+coeff2*wponxj*wponxj
  74.       wcs=wconxj*xj
  75.       arga=wconxj+djonxj
  76.       argc=wponxj/(1.0d0+wponxj)
  77.       argb=dsqrt(1.0d0-argc*argc)
  78.       fshort=1.0d0-xjonxl*(arga*argb-djonxj)
  79.       dwpdvb=xd*dsqdvb
  80.       dadvb=(coeff1+coeff2*(wponxj+wponxj))*dwpdvb*onxj
  81.       dbdvb=-argc*argc*(1.0d0-argc)*dwpdvb/(argb*wps)
  82.       dfsdvb=-xjonxl*(dadvb*argb+arga*dbdvb)
  83.       go to 220
  84. 210   continue
  85.       fshort=1.0d0
  86.       dfsdvb=0.0d0
  87.       wcs=0.05d-6
  88. c
  89. c.....body effect
  90. c
  91. 220   continue
  92.       gammas=gamma*fshort
  93.       fbodys=0.5d0*gammas/(sqphbs+sqphbs)
  94.       fbody=fbodys+fnarrw
  95.       onfbdy=1.0d0/(1.0d0+fbody)
  96.       dfbdvb=-fbodys*dsqdvb/sqphbs+fbodys*dfsdvb/fshort
  97.       qbonco=gammas*sqphbs+fnarrw*phibs
  98.       dqbdvb=gammas*dsqdvb+gamma*dfsdvb*sqphbs-fnarrw
  99. c
  100. c.....static feedback effect
  101. c
  102.       vbix=vbi-eta*vds
  103. c
  104. c.....threshold voltage
  105. c
  106.       vth=vbix+qbonco
  107.       dvtdvd=-eta
  108.       dvtdvb=dqbdvb
  109. c
  110. c.....joint weak inversion and strong inversion
  111. c
  112.       von=vth
  113.       if ( xnfs.eq.0.0d0 ) go to 250
  114.            csonco=charge*xnfs*xl*xw/cox
  115.            cdonco=qbonco/(phibs+phibs)
  116.            xn=1.0d0+csonco+cdonco
  117.            von=vth+vt*xn
  118.            dxndvb=dqbdvb/(phibs+phibs)-qbonco*dsqdvb/(phibs*sqphbs)
  119.            dvodvd=dvtdvd
  120.            dvodvb=dvtdvb+vt*dxndvb
  121.            go to 300
  122. c
  123. c.....cutoff region
  124. c
  125. 250   continue
  126.       if ( vgs.gt.von ) go to 300
  127.       cdrain=0.0d0
  128.       gm=0.0d0
  129.       gds=0.0d0
  130.       gmbs=0.0d0
  131.       if ( icharg.ne.0 ) go to 800
  132.       go to 1000
  133. c
  134. c.....device is on
  135. c
  136. 300   continue
  137.       vgsx=dmax1(vgs,von)
  138. c
  139. c.....mobility modulation by gate voltage
  140. c
  141.       onfg=1.0d0+theta*(vgsx-vth)
  142.       fgate=1.0d0/onfg
  143.       us=uo*fgate
  144.       dfgdvg=-theta*fgate*fgate
  145.       dfgdvd=-dfgdvg*dvtdvd
  146.       dfgdvb=-dfgdvg*dvtdvb
  147. c
  148. c.....saturation voltage
  149. c
  150.       vdsat=(vgsx-vth)*onfbdy
  151.       vpof=vdsat
  152.       if ( vmax.gt.0.0d0 ) go to 310
  153.       dvsdvg=onfbdy
  154.       dvsdvd=-dvsdvg*dvtdvd
  155.       dvsdvb=-dvsdvg*dvtdvb-vdsat*dfbdvb*onfbdy
  156.       go to 400
  157.   310 vdsc=xl*vmax/us
  158.       onvdsc=1.0d0/vdsc
  159.       arga=(vgsx-vth)*onfbdy
  160.       argb=dsqrt(arga*arga+vdsc*vdsc)
  161.       vdsat=arga+vdsc-argb
  162.       dvsdga=(1.0d0-arga/argb)*onfbdy
  163.       dvsdvg=dvsdga-(1.0d0-vdsc/argb)*vdsc*dfgdvg*onfg
  164.       dvsdvd=-dvsdvg*dvtdvd
  165.       dvsdvb=-dvsdvg*dvtdvb-arga*dvsdga*dfbdvb
  166. c
  167. c.....current factors in linear region
  168. c
  169. 400   continue
  170.       vdsx=dmin1(vds,vdsat)
  171.       if ( vdsx.eq.0.0d0 ) go to 900
  172.       cdo=vgsx-vth-0.5d0*(1.0d0+fbody)*vdsx
  173.       dcodvg=1.0d0
  174.       if (vds.lt.vdsat) dcodvd=-dvtdvd-0.5d0*(1.0d0+fbody)
  175.       dcodvb=-dvtdvb-0.5d0*dfbdvb*vdsx
  176. c
  177. c.....normalized drain current
  178. c
  179. 410   continue
  180.       cdnorm=cdo*vdsx
  181.       gm=vdsx
  182.       gds=vgsx-vth-(1.0d0+fbody+dvtdvd)*vdsx
  183.       gmbs=dcodvb*vdsx
  184. c
  185. c.....drain current without velocity saturation effect
  186. c
  187.       cd1=beta*cdnorm
  188.       beta=beta*fgate
  189.       cdrain=beta*cdnorm
  190.       gm=beta*gm+dfgdvg*cd1
  191.       gds=beta*gds+dfgdvd*cd1
  192.       gmbs=beta*gmbs
  193. c
  194. c.....velocity saturation factor
  195. c
  196.       if ( vmax.eq.0.0d0 ) go to 500
  197.       fdrain=1.0d0/(1.0d0+vdsx*onvdsc)
  198.       fd2=fdrain*fdrain
  199.       arga=fd2*vdsx*onvdsc*onfg
  200.       dfddvg=-dfgdvg*arga
  201.       dfddvd=-dfgdvd*arga-fd2*onvdsc
  202.       dfddvb=-dfgdvb*arga
  203. c
  204. c.....drain current
  205. c
  206.       gm=fdrain*gm+dfddvg*cdrain
  207.       gds=fdrain*gds+dfddvd*cdrain
  208.       gmbs=fdrain*gmbs+dfddvb*cdrain
  209.       cdrain=fdrain*cdrain
  210.       beta=beta*fdrain
  211. c
  212. c.....channel length modulation
  213. c
  214. 500   continue
  215.       if ( vds.le.vdsat ) go to 700
  216.       if ( vmax.eq.0.0d0 ) go to 510
  217.       if (alpha.eq.0.0d0) go to 700
  218.       cdsat=cdrain
  219.       gdsat=cdsat*(1.0d0-fdrain)*onvdsc
  220.       gdsat=dmax1(1.0d-12,gdsat)
  221.       gdoncd=gdsat/cdsat
  222.       gdonfd=gdsat/(1.0d0-fdrain)
  223.       gdonfg=gdsat*onfg
  224.       dgdvg=gdoncd*gm-gdonfd*dfddvg+gdonfg*dfgdvg
  225.       dgdvd=gdoncd*gds-gdonfd*dfddvd+gdonfg*dfgdvd
  226.       dgdvb=gdoncd*gmbs-gdonfd*dfddvb+gdonfg*dfgdvb
  227. c
  228.       emax=cdsat*onxl/gdsat
  229.       emoncd=emax/cdsat
  230.       emongd=emax/gdsat
  231.       demdvg=emoncd*gm-emongd*dgdvg
  232.       demdvd=emoncd*gds-emongd*dgdvd
  233.       demdvb=emoncd*gmbs-emongd*dgdvb
  234. c
  235.       arga=0.5d0*emax*alpha
  236.       argc=xkappa*alpha
  237.       argb=dsqrt(arga*arga+argc*(vds-vdsat))
  238.       delxl=argb-arga
  239.       dldvd=argc/(argb+argb)
  240.       dldem=0.5d0*(arga/argb-1.0d0)*alpha
  241.       ddldvg=dldem*demdvg
  242.       ddldvd=dldem*demdvd-dldvd
  243.       ddldvb=dldem*demdvb
  244.       go to 520
  245. 510   continue
  246.       delxl=dsqrt(xkappa*(vds-vdsat)*alpha)
  247.       dldvd=0.5d0*delxl/(vds-vdsat)
  248.       ddldvg=0.0d0
  249.       ddldvd=-dldvd
  250.       ddldvb=0.0d0
  251. c
  252. c.....punch through approximation
  253. c
  254. 520   continue
  255.       if ( delxl.le.(0.5d0*xl) ) go to 600
  256.       wcs2=wcs*wcs
  257.       delxl=xl-(xl**2/(4.0d0*delxl))
  258.       arga=4.0d0*(xl-delxl)**2/xl**2
  259.       ddldvg=ddldvg*arga
  260.       ddldvd=ddldvd*arga
  261.       ddldvb=ddldvb*arga
  262.        dldvd= dldvd*arga
  263. c
  264. c.....saturation region
  265. c
  266. 600   continue
  267.       dlonxl=delxl*onxl
  268.       xlfact=1.0d0/(1.0d0-dlonxl)
  269.       cdrain=cdrain*xlfact
  270.       diddl=cdrain/(xl-delxl)
  271.       gm=gm*xlfact+diddl*ddldvg
  272.       gds0=gds*xlfact+diddl*ddldvd
  273.       gmbs=gmbs*xlfact+diddl*ddldvb
  274.       gm=gm+gds0*dvsdvg
  275.       gmbs=gmbs+gds0*dvsdvb
  276.       gds=gds0*dvsdvd+diddl*dldvd
  277. c
  278. c.....finish strong inversion case
  279. c
  280. 700   continue
  281.       if ( vgs.ge.von ) go to 750
  282. c
  283. c.....weak inversion
  284. c
  285.                 onxn=1.0d0/xn
  286.                 ondvt=onxn/vt
  287.                 wfact=dexp( (vgs-von)*ondvt )
  288.                 cdrain=cdrain*wfact
  289.                 gms=gm*wfact
  290.                 gmw=cdrain*ondvt
  291.                 gm=gmw
  292.                 if (vds.gt.vdsat) gm=gm+gds0*dvsdvg*wfact
  293.                 gds=gds*wfact+(gms-gmw)*dvodvd
  294.                 gmbs=gmbs*wfact+(gms-gmw)*dvodvb
  295.      1                         -gmw*(vgs-von)*onxn*dxndvb
  296. c
  297. c.....charge computation
  298. c
  299.   750 continue
  300.       if (icharg.eq.0) go to 1000
  301.       if (vgs.le.vth) go to 800
  302.       call mqspof(vds,vbs,vgs,vpof,vdsat1,vth,vbin,gamasd,
  303.      1   qg,qc,qb,cggb,cgdb,cgsb,cbgb,cbdb,cbsb)
  304.       go to 2000
  305. c
  306. c.....charge computation for vgs<vth
  307. c
  308. 800   continue
  309.       xqc=xqco
  310.       call mosq3(vds,vbs,vpof,vdsat1,vth,vbin,gamasd,cox,phi,
  311.      1   qg,qc,qb,cggb,cgdb,cgsb,cbgb,cbdb,cbsb)
  312.       qspof=0.0d0
  313.       go to 2000
  314. c
  315. c.....special case of vds=0.0d0
  316. c
  317. 900   continue
  318.       beta=beta*fgate
  319.       cdrain=0.0d0
  320.       gm=0.0d0
  321.       gds=beta*(vgsx-vth)
  322.       gmbs=0.0d0
  323.            if ( (xnfs.ne.0.0d0).and.(vgs.lt.von) )
  324.      1          gds=gds*dexp((vgs-von)/(vt*xn))
  325.       if (icharg.eq.0) go to 1000
  326.       call mosq3(vds,vbs,vpof,vdsat1,vth,vbin,gamasd,cox,phi,
  327.      1   qg,qc,qb,cggb,cgdb,cgsb,cbgb,cbdb,cbsb)
  328. 1000  qspof=0.0d0
  329. c
  330. c.....done
  331. c
  332.  2000 return
  333.       end
  334.